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INTRODUCTION 



I . 

The tracking problem today is of great interest particu- 
larly in the military area. The term "tracking" is considered 
to be the processing of measurements obtained from a target 
in order to maintain an estimate of its current state. The 
latter usually contains: 

a. Kinematic components (position, velocity, acceleration, 
etc . ) . 

b. Other components (radiated strength of signal, spec- 
tral characteristics, etc.). 

c. Constant or s lowly- vary ing parameters (propagation 
velocity, coupling coefficients, etc.). 

A track is a state trajectory estimated from a set of 
measurements which are associated with the same target. 

Measurements are noise-corrupted observations related 
to the state of the target. For example a measurement can 
be direct estimate of position, range, bearing, time differ-’ 
ence of arrival, frequency of narrow band signal, or frequency 
difference of arrival. 

Many approaches have been developed for the solution of 
the tracking problem. Many of them are very reliable and 
work with satisfactory results. The problem however is 
complicated, so there is always need for further improvement 
of the tracking procedures . 
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The complexity of the tracking problem is due to the 
following main reasons: 

a. Presence of countermeasures (evasive maneuvers, 
decoys , etc . ) . 

b. Need of advanced information processing techniques 
for manipulation of the numerous and complex data to 
be possible. 

c. Uncertainty associated with the measurements related 
to the origin of the measurements, i.e., a measure- 
ment may not have originated from the target of interest. 

d. Inaccurate measurements because of random noise. 

e. Random false alarms in the detection process. 

f. Clutter because of reflectors or radiators near the 
target of interest. 

g. Interfering targets. 

The basic components of a tracking system are shown in 
Figure 1.1. 

The term "filter" is often applied to devices or systems 
which process incoming signals or other data in such’ a way 
as to eliminate noise, to smooth signals, or to predict the 
input signal from instant to instant. There is much litera- 
ture covering the theories of estimation, identification, 
modeling, prediction, etc. The design of such filters falls 
in the domain of optimal filtering, which originated with 
the work of Wiener [Ref. 1] and was extended by the work of 
Kalman-Busy [Ref. 2] and others. It was only around 1970 
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that the systematic treatment of tracking multiple targets 
in the presence of false alarms using Kalman filtering 
techniques has started with the work of [Ref. 3] and [Ref. 4] . 

It has been generally accepted that the Kalman filtering 
technique is one of the most desirable of existing tracking 
algorithms. The reason is that it is able to accept many 
inputs (bearing, range, course, speed, doppler, etc.). The 
inputs can be available from a variety of locations and 
sensors. The Kalman filter then weights them properly with 
respect to their errors and generates an estimate of target 
position, course, and speed. It contains also information 
within its structure sufficient to give the user or decision 
maker a useful indicator of the estimate's quality. 

The purpose of this paper is to describe the target 
tracking problem under noisy conditions. The Kalman filter 
is used as the basic tool to treat this problem. Some char- 
acteristic problems are presented and solved. 

In Chapter II some basic concepts and definitions are 
provided from estimation theory. The contents of this chap- 
ter form the theoretical background needed for the rest of 
the chapters. It is assumed that the reader is familiar 
with basic probability theory. 

Chapter III contains a general description of the Kalman 
filter (discrete and continuous) and the nonlinear case 
where a linearization of target tracking problem takes place 
for the application of Kalman filter to be possible (Extended 
Kalman Filter) . 
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In Chapter IV, a brief analysis of the error covariance 
matrix and the information contained in it has been made. 

Also its use in setting error ellipsoids (confidence regions) 
about an estimated target position has been described. 

In Chapter V the maneuvering target problem is described 
and one of the various existing approaches is analyzed, 
named "The White Noise Model with Adjustable Level." The 
same method is used later in Chapter VI, in a relative exam- 
ple. Chapter VI contains three characteristic examples 
covering the application of Kalman filter in some individual 
cases of target tracking problems. Specifically the presented 
examples are the following: 

a. Linear target problem (measurements of target’s posi- 
tion available) . 

b. Non-linear target problem (measurements of target's 
bearing available). 

c. Linear maneuvering target's problem (measurements of 
position available) . The target maneuvers two times 
with different acceleration each time. 

The relative calculations are provided analytically 
together with the corresponding computer programs . The 
resulting plots are also provided with proper comments. 



18 



II. SOME DEFINITIONS AND BASIC CONCEPTS 
FROM ESTIMATION THEORY 



The following definitions and concepts from probability 
and estimation theory respectively are described in detail 
in Reference 5 and Reference 6. 

A. DEFINITIONS 

1 . Bayes 1 Rule for Random Variables 

For random variables Bayes' rule is given by: 



P (x|y) 



P (y I x) P (x) _ P(yx)P(x) 

P (y) J P (y , x) P (x) dx 



( 2 . 1 ) 



The conditional probability of an event is sometimes referred 
to as "posterior" while the unconditional one as "prior. " 

In this case p (x) is called the "prior" probability density 
function (p.d.f.) and p(x|y) is the posterior one. 

2 . Gaussian Random Variables 

The p.d.f. of a Gaussian random variable is 

_ 2 

p(x) = N(x;x,a 2 ) = ~ ^ exp{- — ^ — } (2.2) 

/2tto 2o 



where N( ) represents the normal density with argument x, 
mean x and variance c. 

An equivalent expression is 

x ~ N (x , a 2 ) (2.3) 
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where x is normally distributed with respective mean and 
variance. 

3 . White Noise 

A discrete Gaussian random process v is called "white 
noise" if for any k time points t^ , t^ , . . . , t^_ , the random 
vectors v ( t^ ),..., v ( t^ ) (which are Gaussian) are independent. 
The autocorrelation of this random process is zero for any 
two different times. 

A discrete Gaussian random process is often used as 
a model for measurement noise and disturbance noise. 

4 . Markov Processes 

The common property of Markov processes is the 
following : 

p[x(t)|x(x), x <_ t-jJ = p [x (t) | x ( t^) ] V t > t 1 (2.4) 

i.e., the past up to t, is completely described by the value 
of the process at t. 

The state of a dynamic system with input a white- noise 
x(t) = f[x(t),n(t)] (2.5) 

is a Markov process. 

5 . Random Sequences 

A random sequence, or a discrete time stochastic 
process, is a time-indexed sequence of random variables. 
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x k =’ { x ( j ) , j = 1 , • • • / k } k = l,.... (2.6) 

A random sequence is Markov if 

p(x(k)|x^] = p[x(k)|x(j)] V k > j (2.7) 



The sequence v(j), j = 1,..., is a (discrete time) 
white noise if 



E [ v (k) v ( j ) ] = q5 



kj 



2 . 8 



where 6, . is the Kronecker delta function: 



kj 



1 if k = j 



0 if k j 



(2.9) 



The state of a dynamic system excited by white noise 



x (k+1) = f [k,x(k) ,v(k) ] (2.10) 

is a discrete-time Markov process or Markov sequence. 

A stochastic process x(t) , t e I is a Gaussian white 
(or independent) process if for any m time points t^,...,t m 
in I (m = integer) , the m random n vectors x(t^) , . . . ,x(t m > 
are independent Gaussian random vectors. 

The state of a linear dynamic system excited by white 
Gaussian noise 
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x (k+1) 



( 2 . 11 ) 



= F x (k ) + v (k) 

is a Gauss-Markov process. 

B. BASIC CONCEPTS IN ESTIMATION 

1 • Estimating Problem-Random and Non-random Parameters 
The problem of estimating a (time-invariant) parameter 

x is : 

Given the measurements 

z(j) = h [ j ,x,w( j ) ] , j = 1,... (2.12) 

made in the presence of random noises w(j) , find a function 

~ ~ k 

x (k) = x(k,z ) (2.13) 



where 



z k = {z ( j) , j = 1, . . . ,k} (2.14) 

that estimates the value of x in some sense. 

The function (2.13) is called an estimator and its 
value is called the estimate. We can use two model's for a 
time-invariant parameter: 

a. Non-random: We have an unknown true value x^. Here 
it is desired that the estimates converge in some sense to 
this true value as k -»■ « (non-Bayesian approach) . 
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b. Random: The parameter is a random variable with a 

prior p.d.f. p(x) . In this case a realization of x accord- 
ing to p(x) is assumed to have taken place; It is desired 
for each measurement sequence to yield an estimate that 
converges to the corresponding realization of x and this 
should hold for all x (Bayesian approach) . 

In the second approach given the prior p.d.f. of the 
parameter, its posterior can be obtained using the Bay's 
rule . 

In the first case the likelihood function 

A k (x) = p(Z k ,x) (2.15) 

is used as a measure of how "likely" a parameter value is 
for the observations that are made. 

A usual method of estimation of non-random parameters 
is the maximum likelihood method. Its maximum likelihood 
estimate (MLE) is 



MI, ... . k i . , . . . i 

x (k) = arg max p(Z|x) (2.16) 

x 

The corresponding estimate for a random parameter is 
the maximum a posteriori (MAP) estimate 

''MAP i k k 

x (k) = arg max p (x | Z ) = arg max[p(z jx)p(x)] (2.17) 

x x 
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MAP ML 

The x will coincide with x for a certain prior 

p.d.f., called "diffuse." 

2 . Linear Zstimation--Static Case 

The minimum mean square error (m.m.s.e.) estimate of a 
random variable x in terms of another random variable v 
is the conditional mean E[xjy]. In practice the evaluation 
of the conditional mean is complicated so usually the "linear 
m.m.s.e. estimation" method is applied. According to this 
method the best estimate is such that it is unbiased and the 
estimation error is uncorrelated from the observable ( s ) , 
i.e., they are orthogonal. 

The estimation of two random vectors x and z jointly 
normally (Gaussian) distributed is now examined: 

Let 



Y. = 




( 2 . 13 ) 



The random variable y is normally distributed with mean 



Z 




and covariance matrix (assumed nonsingular) 



( 2 . 19 ) 
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( 2 . 20 ) 



where 



and 




P P 

—xx — xz 



P 



— zx 



p 

— zz 



P. xx = E [ (x - x) (x - x) ' ] 



P = E[(x-x)(z-z')] 

Then the m.m.s.e. estimate of x in terms of z is 



( 2 . 21 ) 



( 2 . 22 ) 



x 



£ 



E [x | z] 



x + P 

— — xz 




( z - z) 



( 2 . 23 ) 



and the corresponding covariance matrix is 



*xx 



/N A 

E [ (x - x) (x - x) 'I z] 



p - p p 1 p 
—xx — xz zz — zx 



( 2 . 24 ) 



We can see that the optimal estimate is a linear function of 
z (this is because of the Gaussian assumption) , while the 
covariance which measures the "quality" of the estimate is 
independent of the observation z_. 

If the random variables are not Gaussian, then the 
"best linear" estimate of x in terms of z (lin. m.m.s.e.) 
is : 
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x = x + P P : (z - z) (2.25) 

— — — xz — zz — 

The m.s.e. matrix corresponding to Equation (2.25) is given 
by 



E [xx 1 ] 



P -P P 1 P = P 
—xx — xz — zz — zx —xx z 



(2.26) 



From the above results it follows that the best estimator 
(in the m.m.s.e. sense) for Gaussian random variables is the 
same as the best linear estimator for arbitrarily distributed 
random variables with the same first and second moments. 

3 . Estimating with the "Least Squares" Method 

The least squares method is another common estimation 
procedure for non-random parameters. For measurements 

z ( j ) = h ( j , x) + w ( j ) (2.27) 

The least square estimate of x is 

k ? 

x LS (k) = arg min J [z(j) - h(j,x)] (2.28) 

x j = l 

Equation (2.28) does not make assumptions about the noises 
w(j). If these noises are independent identical distributed 
(i.i.d.) zero-mean Gaussian random variables, then Equation 
(2.28) coincides with the MLE for these assumptions. 
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The corresponding case for random parameters is the 



minimum mean square error (m.m.s.e.) estimate 



~MMSE ,, , 
a (k) 



^ 2 k 

arg min E[(x -x) | Z' ] 



x 



(2.29) 



The solution to the above is 



xMMSE^j _ e[x|Z^] = / x p(x|Z K ) dx 



\ 2 . 30 



where the expectation is with respect to the conditional 
p . d . f . : 



P(x I Z k ) 



P (Z i x) P (x) 
P (z k ) 



(2. 3i : 



4 . Consistent Estimators 

For a non-random parameter case, the estimator is 
called consistent if the estimate converges to the true 
value in some stochastic sense. Using the convergence in 
mean square, 



lim E { [x (k) - x Q ] 2 } = 0 (2.32) 

k-*<» 

The expectation is 



E [x(Z ) J = 



(2.33) 
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Requirement for the convergence of the estimator 



in the random parameter case (m.s. sense) is 



lim E{[x(k) - x] 2 } = 0 

]^->-co 



The expectation here is over 



and x . 



(2.34) 
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III. KALMAN FILTER AND APPLICATION TO TARGET TRACKING 



A. BASICS ON KALMAN FILTER 

Kalman filtering technique is very popular in target 
tracking applications. It is basically a method of solving 
the problems of optimal filtering and prediction, which 
according to Reference 1 is defined as: 

a. Optimal filtering is the estimation of state vector 
X(t) from data Z_(x) where T <_ t . 

b. Prediction is the estimation of a state vector 'X(t) 
at time t from data Z(x) , where x < t. 

In target tracking both the above problems are encoun- 
tered. The Kalman filter is desirable because as an estima- 
tion model it has the following features: 

a. At time t, the filter generates an unbiased estimate 

A 

X(t) of the state vector X(t) . That means that the 
expected value of the estimate is the value of the 
state vector at time t. 

b. The estimate is a minimum variance estimate (i.e., it 
has the property that its error covariance is less 
than or equal to that of any other linear unbiased 
estimate) . 

c. The filter is recursive and it does not store past 
data . 
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d. The filter is linear which simplifies calculations 



and lends itself to machine computations. Gelb [Ref. 

3] discusses the above features. 

B. DISCRETE KALMAN FILTER 

Assuming that we are dealing with a discrete time system 
of the form: 

X (k+1 ) = o(k+l,k) X (k) + A(k+l,k) u(k) + W(k) (3.1) 

Z(k) = C (k) X (k) + V (k) ' (3.2) 

the following declaration must be done: 

A 

"X( j ! j) " is read as: "X hat of j given j," which means 

that "The estimate of X at time j given measurements at 
times up to and including time j . " 

Then there are eight main components that make up the 
Kalman filter: 

/X 

a. X(k) is the state vector at time k, and X(k) is the 
state vector estimate which is an unbiased minimum 
variance estimate of the true state vector X(k) at 
time k. 

b. The error covariance matrix P(k) is a matrix repre- 
senting the covariance of the difference between the 

/V 

true state vector X(k) and the estimate X(k) , and 
can be expressed as: 
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POO 



(3.3) 



= E { (X(k) - X(k) ) (X (k) _ X(k) ) ’ } 

c. The transition matrix t_(k) is used in the calculation 
of the state vector estimate at the present point, in 
time (k+l|k) , from the state vector estimate of a 

past point in time (kjk) . It is also used in the calcu 
lation of the error covariance matrix at time (k+1 k) , 
using the error covariance matrix at time (k|k) . 

d. The measurement vector Z_(k) , is the data sample or 
observation taken at time k and its components are 
linear combinations of the components of the state 
vector X(k) which have been corrupted by uncorrelated 
Gaussian noise V(k) whose mean is zero and has a 
covariance matrix R(k). 

e. The covariance matrix R(k) associated with the Gaussian 
noise that corrupts the observation or measurement at 
time k, must be provided by the user. It is the covari 
ance of a sensor measurement vector Z_(k) . For example 
a range sensor may have a measurement variance of 60 
yards. A measurement made with that sensor would have 
a scalar value of 60 yards for R(k) . As will be seen 
later, it is important that the value of R ( k ) be 
accurate as it affects the performance of the filter. 

f. The conversion matrix C(k) describes the linear com- 
binations of the components of the state vector X(k) 
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which makes the components of the measurement vector 
Z_(k) . The relationship between the measurement vector 
Z_(k) , the conversion matrix C(k) , and the Gaussian 
noise V(k) is: 

Z(k) = C (k) X (K) + V (k) (3.4) 

g. The Kalman gain matrix G(k) is instrumental in mini- 
mizing the difference between the estimate X(k) and 
the state vector X(k). G(k) is chosen to minimize 
the trace of the estimate error covariance matrix 
P(k) and is used to revise the estimate X(k) of the 
state vector X(k) and the error covariance matrix P(k) 
after the observation Z_(k) has been made. 

h. The last component to be discussed is the system 
dynamics covariance matrix Q(k). The assumption of 
the Kalman model is that the state vector X(t) exists 
in a random environment that is Gaussian with zero 
mean and covariance Q(k) . The system dynamics covari- 
ance matrix must be input by the user, and its 
determination is important as it also affects the 
performance of the filter. 

The remainder of this subsection shows the relation- 
ships between the algorithm. Two important points follow: 
a. The transition matrix g (k) , measurement noise covari- 
ance matrix R(k), systems dynamics covariance matrix 
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Q(k) and the conversion matrix C(k) may be different 
for each point in time k. 

b. The filter must be initialized by the user providing 

/\ 

the initial estimate X(k) and its associated estimate 
error covariance matrix P(k) . A poor initialization 
will require more observations for the algorithm esti- 
mate to converge near the value of the state vector. 

A Kalman filter iteration can be divided into phases: 
Prediction and Filtering. 

During the prediction phase, the state vector esti- 
mate X(k+ljk) and its error covariance matrix P_(k+1 k) are 
updated from the previous value at time (k|k) to the time 
(k+l|k) when the current measurement Z_(k) is observed. 

During the same phase the system dynamics is introduced into 
the algorithm. The covariance matrix Q(k) accounts for the 
system dynamics (environment) from time (k) to (k+l;k). 

/*S 

The update is performed by multiplying the estimate X(k|k) 
and the error covariance matrix P(k|k) by the transition 
matrix £(k+l,k) as follows: 

X (k+1 | k) = £(k+l,k) X (k | k) (3.5) 

P (k+1 | k) = d(k+l,k) P (k | k) d' (k+l,k) + Q(k) (3.6) 

If the system has an input u, the first of the above 
equations is extended as follows: 
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X (k+1 , k) = 2(k+l,k) X (k k) + (k+l,k) u(k) 



(3.7) 



During the filtering phase, the estimate X(k k) of 



the state vector X(k k) and the error covariance matrix 
P(k'k) are revised, based onthe latest measurement Z_(k) 
observed at time k. To do this, the Kalman gain matrix 
G(k) is first computed using the error covariance matrix 
P(k k-1) which was updated in the prediction phase. The 
covariance matrix R(k) which is the covariance of the Gaussian 
noise associated with the latest measurement Z_(k) and the 
conversion matrix C(k) are also associated with the latest 
measurement. The sequence of computations during the filter- 
ing phase is : 



G (k) = P(k|k-1) C' (k) [C(k) P (k|k-l)C’ (k) +R(k)] 



-1 



(3.3) 



X (k | k) = X (k | k-1) + G(k) [Z (k) - C(k) X(klk-l)] 



(3.9) 



P(kik) = [I - G (k) C (k) ] P { k | k-1) 



(3.10) 



The measurement prediction covariance is 



S ( k+1 ) = C(k+l)P(k+l|k)C' (k+1) + R (k+1) 



(3.11) 



The innovation or measurement residual is 
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V (k+1) 



z ( k+1) 



z (k+1 | k) 



(3.12) 



An important property of the innovations is that 
they are an orthogonal sequence, i.e., 

E ^j> = S k « kj (3.13) 

where 6, ■ is the Kronecker delta function. 

K J 

It must be noticed that, at every stage k, the entire 
past is summarized by the "sufficient statistic " x(k k) and 
the associated covariance. The covariance equations are 
independent of the measurements. The covariance equations 
can thus be iterated forward offline. 

Figure 3.1 shows a schematic of the Kalman filter 
algorithm and relationships among the components described 
above. Reference 7 provides a detailed discussion of Kalman 
filtering and optimal estimation. 

C. CONTINUOUS KALMAN FILTER 

In order to go from a discrete to a continuous system 
of the form: 

X(t) = F ( t ) X(t) + B ( t ) W(t) (3.14) 

Z (t) = C (t) X (t) + V(t) (3.15) 
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PEAL WORLD 
(not observable 
by user) 



FILTER 

(Input from Time k) 




Figure 3.1 Schematic of the K.F. Algorithm 



where W,V are zero mean, uncorrelated, white noise processes 
with covariance matrices Q and R respectively, we have to 
consider the situation at the limit as t(k) - t(k-l) = At -*• 0 . 

At this limit we have the following equivalences: 

<p (k) -*• I + F t (3.16) 

Q (k) - B Q B 1 At (3.17) 

R ( k ) - R/At (3.18) 

R(k) = E(V(k) V' (k) ) is a covariance matrix, while: 

R ( t ) = E ( V ( t ) V ' ( T ) ) = R ( t ) 6 ( t - T ) (3.19) 

is a spectral density matrix. The covariance matrix R(t)o(t -t) 
has infinite valued elements. To approximate the continuous 
white noise process by the discrete white noise sequence, 
the pulse lengths (At) may be shrinked and their amplitudes 
may be increased, such that R(k) -*■ R/At. In other words, 
as At -* 0 , the discrete noise sequence tends to one of an 
infinite valued pulse of zero duration, such that the area 
under the "impulse" autocorrelation function is R(k)At, 
equal to the area R under the continuous white noise impulse 
autocorrelation function. 

With the above expressions in mind, the procedure is to 
write the proper difference equations and to observe their 
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behavior in the limit as it - 0. It is assumed that R is 
non-singular, i.e., R 1 exists. 

Finally the continuous Kalman filter equations are sum- 
marized as follows: 

a. System Model 

X(t) = F ( t ) X ( t ) + B ( t ) W(t), W(t) - N(0,Q(t)) (3.20 

b. Measurement Model 

Z(t) = C(t) X (t) + V(t), V(t) ~ N(0,R(t)) (3.21) 

c. Initial Conditions 

/\ A A 

E [ X ( 0 ) ] = X Q , E [ (X ( 0 ) -X Q )(X(0) - X 0 ) ’ ] = Pq (3.22) 

d. Assumptions 

R ^"(t) exists 

e. State Estimate 

• /N A A A 

X ( t ) = F(t) X(t) + G(t) [Z (t) -C(t)X(t) ] , X(0) = X Q (3.23) 

f. Error Covariance Propagation 

P(t) = F(t)P(t) +'P(t)F'(t) + B(t)Q(t)B' (t) -G(t)R(t)G' (t) 

where P(0) = Pq (3.2 4)- 
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g. Kalman Gain Matrix 

G ( t ) = P (t) C' (t) R -1 (t) when E[W(t) V'(i)] = 0 (3.25) 

G(t) = (P(t)C’(t) + B (t) J (t) ] R _1 (t) (3.26) 

For more details, See Reference?, pp . 119-124. 

D. EXTENDED KALMAN FILTER 

In target tracking it is usually needed to estimate the 
present target position, course and speed and to predict 
future target position based on the present estimate. The 
Kalman filter works very well for target tracking since the 
algorithm can provide an unbiased, minimum variance estimate 
of the target's state based on varied observations (filter- 
ing), predict future position using the prediction phase of 
the filter and provide an indicator of the estimate error 
through the estimate error covariance matrix. 

Usually in practice, the state and/or measurement equa- 
tions are not linear. Since the Kalman filter is applied to 
linear cases, it is necessary to find a "method" to use it 
in nonlinear estimation problems. One approach is to derive 
an optimal filter for the nonlinear problem. Another approach 
(more usual) is to linearize the problem and apply Kalman 
filter theory to the linearized problem. The highlights of 
the second method are presented in the following: 
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The system and measurement equations are assumed zo be 
of the form: 

Discrete 



x(k+l) = a (x ( k) , u ( k) , k ) + w(k) (3.27) 

z(k) = c (x (k) ) + v (k) (3.28) 

Continuous 



x(t) = f(x(t),u(t),t)+w(t) (3.29) 

z(t) = h(x(t)) + v ( t ) (3.30) 

Here it is assumed that we deal with a discrete model. 

It is assumed that we have available a trajectory 

( 0 ) 

x (k) , k = 0,1,2,.... The vector function a (x (k) , u (k) , k) 
is expanded in a Taylor series about the nominal trajectory 
x^ (k) . Then the linearized state equations can be written 
as : 

x(k+l) = A(k) x(k) + a(x^) (k) ,u(k) ,k) - A(k) x^^ (k) + w(k) 

(3.31) 

where A(k) is defined to be the first partial derivatives 
of the nonlinear function a(x^^ (k) ,u(k) ,k) , with respect 
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to the state vector x(k), i.e., 



A (k) 



3 a > 

3 -| (x (0) (k) , u (k) , k) 



( 3 . 32 ) 



The accuracy of this approximation depends on the differ- 
ence between x^ (k) and the actual trajectory x(k) . The 
middle two terms are treated as deterministic inputs. 

Now the measurement equation has been considered. We 

have 



z(k) = c (x (k) ) + V (k) 



( 3 . 33 ) 



The nonlinear vector c is expanded again about the nominal 
( 0 ) 

trajectory x (k) . Then the measurement equation can be 
written as 



z(k) = H(k)x(k) + c(x^(k)) - H(k)x^°^(k) + v(k) 



( 3 . 34 ) 



where H(k) is defined as the first partial derivatives of 
the measurement function c(x^^ (k) ) with respect to the 
state vector X(k), i.e., 



A 

H(k) = 33 ? 



x (0) (k) 



( 3 . 35 ) 



As in the linearized state equation, the two middle terms 
are known time-varying quantities which can be handled as if 
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they were a time-varying measurement bias, 
tion it is defined 



For simplifica- 



u’ (k) = a(x (0) (k) ,u(k) ,k) -A(k)x (0) (k) (3.36) 

z' (k) = z(k) - c(x (0) (k)) + H(k)x (0) (k)) 

= z ( k ) - b ( k ) ; 3 . 3 7 , 

so that 

x(k+l) = A(k)x(k) + u' (k) + w(k) (3.38) 

z’ (k) = H (k) x (k) + v(k) (3.39) 

With these linearized equations, the appropriate Kalman 
filter equations are: 

a. The gain equation 

G(k) = P (k | k-1 ) H' (k) [H(k) P (k |k-l) H’ (k) +R(k)] _1 (3.40) 

b. The covariance of estimation error equation 

P (k | k-1 ) = A (k-1) P (k-1 | k-1) A’ (k-1) +Q(k-1) (3.41) 

P (k | k) + [!_ - G (k) H (k) ] P (k | k-1) (3.42) 
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Filter update equation 



c . 



/X /N /v 

X (k | k) = X(kjk-l) + G(k) (z (k) - c (x (k; k-1) ) ] (3.43) 

d. Prediction equation 

A A 

x(k+l|k) = a (x (k j k) , u (k) , k) (3.44) 

A block diagram of the system and estimator is shown in 
Figure 3.2. 

The gains can be pre-computed and stored if it is assumed 
that the nominal trajectory is known in advance. 

To answer the question, "where does the nominal trajectory 
x(k) come from?", two approaches are usually used. 

In the first, an approximate trajectory is used. This 
trajectory is known in advance from known data or may have 
been computer by solving the state equations 

x (0) (k+1) = a (x (0) (k) ,u(k) ,k) (3.45) 

with the initial condition x^ (0) = E[x(0)]. If the uncer- 
tainty in X(0) is large, the solution of the above equation 
may be "far" from x(k) to make the linearization sufficiently 
accurate. If an accurate nominal trajectory is available, 
the gain can be computed off-line and stored. 

In the second approach, the estimates produced by the 
filter are used as the nominal trajectory about which the 
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v(k) 




>DDEL 




Figure 3.2 E.K.F. Block Diagram of the System and 
Estimator 
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linearization is performed. The matrices A(k) and H(k) must 
be used to generate G(k) . The best trajectory information 
available when H(k) must be evaluated is x(kjk-l); when A(k) 

A 

is to be evaluated, x(k'k) is available. 

The filter that results from using this approach is 
called an Extended Kalman Filter. The gain and covariance 
equations must be solved on-line. The procedure of the 
filter's calculations is given below: 

A 

1. Start with x(Oj-l) and evaluate H(0) using: 

3c 

H(k) - < (3.46) 

- (x(k k-1) ) 



2 . 



P ( 0 | - 1 ) = M = E { (x ( 0 | - 1 ) - x (0) ) (x(0 | -1) - x ( 0 ) ) ' } (3.47) 

Use this matrix to compute G(0) given by 
G(k) = P(k|k-1)H' (k) [H (k) P (k I k-1 ) H 1 (k) +R(k)] _1 (3.48) 



3. Compute x(0|0) from: 



X(k|k) = x (k | k-1) + G(k) [z(k) -c(x(k|k-l)] (3.49) 



or 



45 



X ( 0 i 0 ) 



c(x(0 -1) ) 



(3.50) 



= x ( 0 -1) + G ( 0 ) [z(0) - 

and x(l 0) from: 



X (k+1 i k) = a (x(k | k) ,u(k) ,k) 
4. Compute P(0 0) from: 



(3.51) 



P (k i k) + [I - G ( k ) H (k) ] P(k k-1 



3 . 52 



5. Since x(0,0) is available before A(0) is calculated, 

/s. 

the value of x(0 0) is used, hence: 



A ( 0 ) = 



3a 

oX 



(3.53) 



(x ( 0 I 0) , u ( 0 ) ,0) 



and then P(1|0) is computed from: 



P (k k-1) = A(k-l) P (k-1 Ik-1) A' (k-1) + Q(k-l) 



(3.54 



The process is repeated by recycling through steps 1 to 5 . 
The Extended Kalman Filter presented was obtained by a 
Taylor series expansion up to the first order terms. This 
filter obviously introduces errors in the equations where the 
higher order terms are neglected. 

There are several ways of compensating for these errors. 
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Very roughly , an addition of "artificial process 
noise" covariance can be made, for compensation of the errors 
is the state prediction. Also a multiplication of the state 
covariance by a scalar 1 > 1 at every sampling time can be 
carried out. This multiplication is equivalent to the filter 
having a "fading memory," i.e., at every sampling time the 
past data is "discounted" by attaching to it a higher covari- 
ance (lower accuracy) . 

The above subject is covered with more details in Reference 
8, pp. 4-52 — 4.58 and in Reference 5, pp . 3 . 3-1--3 . 4-1 . 
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IV. 



ERROR COVARIANCE MATRIX AND TARGET 
TRACKING QUALITY 



A. DEFINITION OF ERROR COVARIANCE MATRIX 

The error x in the estimate of a state vector is defined 

/\ 

to be the difference between the estimated (x) and -he 
actual (x) values: 



x = x - x (4.1) 

The above difference is known as the error vector or 
estimate error. The covariance of the estimate error is 

P = E ( (X ( t ) - X ( t ) ( X ( t ) — X ( t ) ) ' ] (4.2) 

It provides a statistical measure of the uncertainty in 
x. It is possible to discuss the properties of the covari- 
ance matrix independently of the mean value of the state. 

The information contained is sufficient to generate indicators 
of the estimate quality. 

B. INFORMATION CONTAINED IN THE ERROR COVARIANCE MATRIX 
There are five important characteristics of the error 
covariance matrix which relate the matrix to the state 
vector and its estimate. 

(1) The error covariance matrix of an n-component state 
vector is an n by n symmetric matrix. 
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(2) The diagonal elements of the error covariance matrix 
are the mean square errors of the -error vector 
components . 

(3) The trace of the error covariance matrix is the 
mean square length of the error vector. 

(4) The off-diagonal terms of the matrix are correlations 

A 

between the elements of the error vector X(t) - X(t). 

(5) The error covariance matrix P(k) tends to the system 
dynamics covariance matrix Q(k) , as k goes to 
infinity. This means that as more information is 

'available about the state vector (observations) the 
estimate uncertainty approaches the uncertainty of the 
environment in which the state vector exists. 

C. ERROR ELLIPSOIDS ABOUT A TARGET POSITION 

Frequently it is significant to have available an indica- 
tion of the quality of the estimates. One approach to 
achieve this is the proper use of the error covariance 
matrix P(k|k). The outline of this approach is described 
below. 

It is assumed that the initial state of the model x(0) 
and the random processes v(k), w(k) are Gaussian. If this 
assumption is satisfied then it can be shown that: 

x(k),x(k|k) and e(k|k) = x(kjk) - x(k) 
are Gaussian. The results obtained apply only to this case. 
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The probability density function for e(k k) can be 



written as 



f E (e(k : k) ] 




(4.3) 



For a fixed time k this expression can be written as 




1 




(2tt) n/2 P 1/ 2 



1 -1 

where w = p (k.k) . In order to determine the locus of 



points where the density function f_(e) , has a constant 

h — 

value the above equation has to be examined. It is seen 
that fjr(e) has a constant value whenever 



It can be shown that the points e which satisfy Equation 
(4.5) are hyperellipsoids (in two dimensions, ellipses). 

If the left side of (4.5) is expanded for the two-dimensional 
case (the same approach can be extended to n dimensions), 
we have: 



e 1 we 



constant 




(4.5) 




2 

e . + w 



1 

2 




1 



2 



(4.6) 
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where the symmetry of w has been used (w^ 2 = w 21^ ' Equation 

2 

(4.6) is an ellipse (w^ > 0, w 22 > 0 and w 22 > w 12 ). 

It is not easy to be recognized as such because its principal 
axes do not coincide with the coordinate axes as shown in 
Figure 4.1. 

First the principal axes must be determined and then the 
ellipse can be rewritten in terms of and y/ ^ as 

coordinate axes . 




Figure 4.1 Error Ellipse 



The ellipse in the new coordinate system is described by 



,2 ,2 

X l y i + A 2 y 2 



y- 



l/A 



1/a. 



= 1 ' 



(4 . 7) 



where y^ and y^ are eigenvectors of w and ^ and ^ are 
the corresponding eigenvalues. 

As it has already been mentioned, W = P . Equation ( 4.“/ 
is an ellipse in terms of the eigenvalues ^ and 2 and the 
eigenvectors y^' 1 and Q f w. The expression which gives 

the ellipse in terms of the eigenvalues and eigenvectors 
of P is 




where y/ ^ is an eigenvector for P and a^ = 1 /a^ is the 
corresponding eigenvalue. 

This result can be generalized to n dimensions as well. 
Given the quadratic form 



1 

2 



e 1 




(4.9) 



the eigenvalues of P are a^a^; . . . ,a and the corresponding 
eigenvectors are y ,y , . . . . The quadratic form 
(Equation (4.9)) describes a hyperellipsoid which can be 
written as 
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(4.10) 




All vectors e in the n-dimensional space can be written as a 

linear combination of the eigenvectors y ,y , . . . ,y_ 

The coefficient of in the linear combination is y^ , 

( 2 ) 

the coefficient of y_ is y^, etc. 

So the surfaces of equal probability density are 
ellipses (or hyperellipsoids). 

The problem can now be stated as follows: 

For a specified value of A, what is the probability that 
e lies within or on the ellipsoid? 

These probabilities have been tabulated below for a 
few values of n and A. 

TABLE I 

Probabilities for Error Ellipsoids 

A 

n 1 2 3 



1 


.683 


. 955 


. 997 


2 


. 394 


.865 


. 989 


3 


.200 


. 739 


. 977 . 


For example, 


for a 


system having 


n = 2: 



Probability error inside A = 1 ellipse = 0.394 
Probability error inside A = 2 ellipse = 0.865 
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Probability error inside A 



3 ellipse = 0.989 



Therefore, if the error covariance matrix P is given, the 
error ellipsoids can be determined by finding the eigenvalues 
and eigenvectors. 

The error ellipses are useful in visualizing the estima- 
tion error. By using them we can consider the true state 
value to lie within a certain region surrounding the estimate. 

For details and derivations, see Reference 8, pp . 4.44-- 

4.51. 
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V. MANEUVERING TARGETS 



A. GENERAL DESCRIPTION OF THE PROBLEM 

The maneuvering target is generally described by the 
equation : 

x (k+1) = F (k) x (k) + G (k) u (k) + v(k) (5.1) 

where the matrices F(k) , G(k) are assumed known and the 
process noise v(k) is zero mean, white random sequence with 
covariance matrix Q(k) . The main characteristic of the 
maneuvering target equation is that the inputs u(k) are 
unknown . 

In the following, linear models are examined for sim- 
plicity but the same techniques can be applied to nonlinear 
cases . 

A number of different approaches to the maneuvering target 
problem have appeared’ in the literature. The most commonly 
used model is, due to simplicity requirements, a kinematic 
model with states consisting of position, velocity, and in 
most cases also acceleration, driven by white noise. 

It is possible to divide the different approaches into 
two broad categories: 

A. The unknown input (maneuver command) is modeled as a 
random process . 

B. The unknown input is estimated in real time. 
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The random process type models can be further classified 
into two categories, depending on their statistical properties 
Al . White noise 

A2 . Autocorrelated (Markov) noise 

All these methods are approximations because in general, 
the maneuvers are not stochastic processes. 

In the input estimation case the assumption of having a 
constant input over a certain period of time is usually made. 
The estimation can be based on the least squares criterion 
and can be used in the following ways: 

31. The estimated input corrects the state estimate. 

B2 . The input becomes an extra state component (state 

is augmented) that is reestimated within the state. 

In the following section only one method is described. 

This method is illustrated by an example in Chapter VI. It 
was selected from the others, due to its simplicity. 

For detailed descriptions of different approaches to the 
problem, see Reference 5. 

B. THE WHITE NOISE MODEL WITH ADJUSTABLE LEVEL 

In this method a certain low level of process noise is 
assumed in the filter. 

A maneuver is considered as a large innovation. The 
detection procedure is based on the square norm of the 
innovations 



e v (k) = V' (k) S 1 (k) V (k) 



(5.2) 
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where 



V(k) = z (k) - z(k!k-l) 



(5.3) 



Based on the target model (for a non-maneuvering target) 
a threshold is set up 



where a is arbitrary. For example, a = 0.05. 

If the threshold is exceeded, Q(k-l) is multiplied by 
a scaling factor <p until e^. is reduced to the threshold 



When using the factor the covariance of the innovations 
takes the form: 

S (k) = C(k) (F(k-l)P(k-l | k-1 ) F’ (k-1) + <j)Q(k-l) ] C’ (k) + R(k) 

(5.5) 

This obviously reduces the value of e (k) . 

An analogous technique is possible to be used to lower 
the process noise level after the maneuver. 




max 



1 - a 



(5.4) 



C 



max 



57 



VI. CHARACTERISTIC COMPUTER EXAMPLES 



A. A LINEAR KALMAN FILTER EXAMPLE 
1 . Problem Description 

The target is assumed to be described by the system: 



x ( k+1 ) = 



*1 


- 

T 


' i 

x (k) + 


T 


0 


1 




1 

L J 



w (k) 



k = 0,1, 



( 6.1 



The available measurements are of the type: 



z (k) = [1 0 ] x (k) + v ( k) k = 1,2, .. . 



( 6 . 2 ) 



where w(k) ~ N(0,Q), v(k) ~ N(0,R) are mutually independent, 
zero mean, white random sequences. 

The initial state is 



x ( 0 ) 



5 

1 



(6.3) 



Two measurements z(0) and z(l) are made to initialize 
a Kalman filter. 

The sampling time is given as T = 0.1. 

The process and measurement noise are to be yielded 
by a Gaussian random number generator, given that Q = R = 0.02. 
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Case 1 



After the run of the Kalman filter for k = 2 ,..., 100, 
the following expressions are useful to be plotted: 

a. True trajectory vs. estimated trajectory in the x^,x 2 
plane . 

b. Position error variance 

PudlD , P L ]_ (2 i 1 ) , P 11 (2|2) . 

c. Normalized position error 

[x L (k) -x 1 (k k) ] / [p L1 (k k)] 1/2 , k - 2,. ..,100 

d. Velocity error variance 

p 22 (l|l), P 2 2 (2 | 1 ) / p 2 2 ( 2 | 2 ) , . - . 

e. Normalized velocity error 
[x 2 (k) - x 2 (k | k) ] / [p 22 (k | k) ] 1//2 

f. Normalized state error squared 

[x (k) - x (k | k) ] ' P 1 (k) [x (k) -x(k|k)] 

g. Normalized innovation error 

[z(k) - z (k|k-l) ]/[ Pll (k|k-l) + R] 1/2 
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Case 2 



In order to see the effect of running the Kalman 
filter with a different Q than that of the model, it is helpful 

to plot the same expressions as in Case 1, using a different 

-2 

Q for the target's model, say Q = 9 * 10 
Case 3 

Finally, it is also interesting to see the change 
on the same expressions as the Q of both the model and the 

filter increases. For this reason new plots of these expres- 

-2 

sions are to be obtained, using Q = 9 x 10 for both the model 
and the filter. 

2 . Analysis 

a. True Trajectory 

The matrix form of the given model is converted 
in the following equations: 



x^ (k+1) 



x-^ ( k ) + Tx 2 ( k ) + Tw ( k ) 



(6.4) 



x 2 (k+1) 



x 2 ( k ) + w ( k ) 



(6.5) 



where T = 0.1, w(k) is generated by a random number generator 



function and 



5 



x ( 0 ) 



1 
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b. Measurements 

From the given measurement equation it is 

obvious that 



z_(k) = x 1 (k) + v (k) 



(6.6) 



where x^(0) = 5 and v(k) is generated by a random number 
generator function. 

c. Derivation of x(ll), p ( 1 1 ) 

It is reasonable to start with the following 
estimations for the position and velocity vectors: 



x(l 1) 



then 



x 1 (lj 1) 

x 2 (1 | 1) 



z (1) 

z (1) - 2 ( 0 ) 
T 



(6.7) 



x 1 (l|l) 



x 2 ( 1 | 1) 



x x (l ) + v ( 1 ) 

X 1 (1) - x x (0) 

T 



v (1) - v ( 0 ) 

T 



(6.8) 



(6.9) 



matrix 



Using these values in the initial error covariance 



P ( 1 | 1) 



X x (l|l) x 1 (l|l)x 2 (l|l) 

x 1 (l|l)i 2 ( 1 1 1 ) x 2 (l|l) 



( 6 . 10 ) 
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where 



x 1 (l|l) = x^l 1 ) 

X 2 (1 j 1) = X 2 (1 | 1) 

It is obtained that 



x 1 (l 1) 
x 2 (ljl) 



p (1 i 1) 




R 

T 



2R 

T 



( 6 . 11 ) 



( 6 . 12 ) 



(6.13) 



d. Run of the Kalman Filter 

(1) Calculation of Gains . The known equation 
of the Kalman Filter is: 



P ( k | k-1 ) = b P (k-1 I k-1 ) y + Q' (k-1) 



(6.14) 



It has been converted in matrix form and after some manipula- 
tion, the following equations are obtained: 



P 11 (k|k-1) = P ] 2 (k-1 | k-1 ) + TP 21 (k-1 |k-l) 



-11 



21 



+ TP, - (k-1 | k-1) + T 2 P 22 (k-1 | k-1 ) + T 2 Q(k)) 



(6.15) 



P 12 (k|k-1) = P n 0 (k-1 | k-1 ) + TP,, (k-1 | k-1) + TQ (k) (6.16) 



-12 



-2 2 



P 21 (k|k-1) = P ? 1 (k-1 | k-1) + TP 22 (k-1 | k-1 ) + TQ(k) (6.17) 



-21 



2 2 
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—2 2 ( k j k — 1 ) = P_ 22 (k-1 1 k-1) + Q(k) (6.18) 

Keeping the same procedure, the Kalman Filter 

equation 

G (k) = P(k|k-l)C' [C P (k|k-l)C' + R(k) ] _1 ' (6.19) 

takes the following form 

P u (k j k-1) 

G l (k) = P 11 (k|k-1) + R ( k) (6.20) 

Finally the equation 

P (k | k) = [I - G(k)C]P(k|k-l) (6.22) 

is analyzed in the following: 

P 1]L (k|k) = P i;L (k|k-l) -G 1 (k)P 11 (k|k-l) (6.23) 

P 12 ( k|k) = P 12 (k|k-l) -G 1 (k)P 12 (k|k-l) (6.24) 

P 21 (k|k) = P 21 Ck | k-1) -G 2 (k)P 11 (k|k-l) (6.25) 

— 2 2 ( k | k ) = P 22 (k|k-1) - G 2 (k) P (k | k-1 ) (6.26) 



63 



2 . Estimated Trajectories 



Following the procedure described above the estimated 
trajectories given by the equations: 



x(k,k-l) = x (k-1 k-1) + A u (k) 


(6.27) 


/s ^ ^ 

x(kk) = x(kjk-l) + G(k) [z (k) - C(k)x(k k-1)] 

take the form: 


(6.28) 


x^(kjk-l) = x^(k-l|k-l) + Tx^ (k-1 k-1) 


(6.29) 


i — 1 

i — 1 

1 

CM 

< X| 

II 

i — l 

l 

CM 

< XI 


(6.30) 


A A 

x x ( k ] k ) = G 1 (k)£(k) + (l-G(k) ) x (k | k-1) 


(6.31) 


A A /N 

x 2 ( k | k ) = x 2 (k ( k-1 ) + G 2 (k)z(k) - G 2 (k) x x (k | k-1 ) 


(6.32) 


e. Normalized State Error Squared 

The normalized state error squared is given 
the expression: 

/V -1 /■s 

[x(k) -x(k|k)]'P (k) [x(k) -x(k|k)] 


by 



By manipulating this expression in its matrix 
form it is obtained that it is equal to the following: 
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P 11 P 22 “ P 21 P 12 t t x i-(’<) -Xi (k | k) ] P 22 + [x >2 (k) -x 2 (k j k) ] P.^ 

A A 

- 2P 21 [x 2 (k) x 1 (k) - x 2 (k) (k j k) -x 2 (k k)x^(k) 

/V A 

+ x 2 (k|k)x 1 (k|k) ] ] 

f. Normalized Innovation Error 
By making use of the equation 

/V # /N 

£(k|k-l) = cx(k|k-l) (6.33) 



the given expression for the innovation error becomes 



[z (k) - z (k 1 k-1 ) ] 
[P 11 (k|k-1) + R] 1/2 



[z (k) - x x (k | k-1 ) ] 
[P u (k|k-1) + R] 1/2 



(6.34) 



All the above have been properly set in a FORTRAN 
program (pp. 178~180)- The results corresponding to cases 
1,2, and 3 are indicated on pages 66-89, respectively. 

Some comments have been made on the computer program outputs 
on pages 65 and 90-94. 

3 . Comments on the Graphs 

a. Normalized Expressions c, e, g of Cases 1, 2, and 3 
The normal distribution has p.d.f.: 

1 , —.2 2 
1 (x - x) /a 

f (x) - — — e (6.35) 

a/ 2tt 
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Figure 6.2 Estimated Trajectory in 
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Figure 6.3 Position Error Variance. Case 1(b) 
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Figure 6.4 Normalized Position Error. 
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Figure 6.5 Velocity Urror Varianc 
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Figure 6.6 Normalized Velocity Error. 
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Figure (>.7 Normalized 
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Figure 6.8 Normalized Innovation Error. 
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Figure 6.9 True Trajectory in 
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Figure 6.10 Estimated Trajectory in x.,x Plane. Case 2(a) 
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Figure 6.11 Position hlrror Vnrinnco. tuso 2(h) 
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Figure 6.12 Normalized Position Error. 
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Figure 6.13 Velocity Error Vai iunc<>. (Vise* 2(d) 
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Figure 6.14 Normalized Velocity lirror. 
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Normal lzed St a to Krror 



— gj«- 




81 



Figure 6.16 Normalized Innovation 
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Figure 6.17 
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Figure 6.18 Estimated Trajectory in x,,x I'Lane. 
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Figure 6.20 Normalized Position Urror. Case 
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l*' i cjure (>.21 Velocity 
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Figure 6.22 Normalized Velocity Error. Case 3(e) 
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Figure 6.23 Normalized State Krror Squared. 
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Figure 6.24 Normalized Innovation Frror. Case 3(g) 



By standardize the normal distribution curve we 
can have a single curve that may be adapted to all values of 
the mean as well as differing values of the standard deviation. 
This standardization may be accomplished by substituting 
z = (x-x)/o in Equation (6.35). 

The above concept is the case for the subpara- 
graphs c, e, g. In Case c, for example, 

z = [x 1 (k) -x 1 (k,k)]/[P 11 (k;k)] 1/2 . 

The standard normal distribution curve has a 
mean of zero and variance and standard deviation equal to 
one . 

By integrating Equation (6.35) (after the substi- 
tution of z) , it is found that: 

Between -la < z < +la, 68.2% of the area under the curve 
is included. 

Between -2a < z < +2a, 94.5% of the area under the curve 
is included. 

Between -3a < z < +3a, 99.74% of the area under the curve 
is included. 

With the above in mind, the probability that the 
normalized position, velocity and innovation errors lie between 
±2a is 94.5% (or between ±3a is 99.74%) . 

From the graphs it is seen that in the two 
cases where the value of Q is common for the model and the 
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filter, the three normalized errors lie 100% between 
z 3c . 

In the case where Q is different between the model 
and the filter it is seen that these errors exceed at some 
points ±3o. But generally it seems that the law of 99.74% 
between ±3c is applied satisfactory. 

b. True Trajectory vs. Estimated Trajectory 

(expr. a) 

( 1 ) Q = 0.02 for Model and Filter (pp. 66-67) . 

There is considerable difference between true and estimated 
curve at the beginning but as time progresses, there is a 
continuous improvement. Near the end, the two curves nearly 
coincide . 

( 2 ) Q = 0.02 for the Filter and Q = 0.09 for 
the Model (pp. 74- 75) . There is a significant difference 
between true and estimated curve at the beginning which is 
slightly improved at the end of the curves. 

( 3 ) Q = 0.09 for the Model and Filter (pp. 82-83) . 
Same as in (1) , with the only difference that at the end of 

the curves, the improvement is not nearly a coincidence but 
it is better than that of Case (2) . 

Generally the rate of the improvement is 
greater' at small values of time because the gain is inversely 
proportional to time (G(k) = 1/k+l) and it weights the correc- 

A 

tion term [£(k) - C(k)x(k|k-1) less heavily as time progresses. 

The fact that the improvement between true 
and estimated curves is less in Case (3) than that of Case (1), 
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was expected. The reason is than Q in Case ( 3 ) has a greater 
value and as it is known, Q increases the uncertainty in che 
one seep prediction: 

?(k k-1) = P (k-1 1 k-1) eg (6.36) 

This affects also the P (kf k) : 

P (k k) = [I - G(k) c (k) ]P (k k-1) .6.3" 

The fact that the estimated curve in Case 

(2) differs from the true curve more than that of Cases (1) 

and (3) and also that the improvement is not so significant 

as in the other cases was also expected. The reason is than 

the Q's of the mode-and the filter are different, so this 

difference affects the P and G of the filter negatively, 

resulting in a difficulty in approaching the true trajectories. 

c. Normalized State Error Squared Error (expr. f) 

It is known from the theory that (x -x) 'P ^ (x - 

is the sum of the squares of n independent zero-mean, unity 

variance Gaussian random variable, i.e., N(0,1). This means 

that (x-x)'P 1 (x - x) has a chi-square distribution with n 

degrees of freedom (n is the dimension of vector x) . 

2 

If the chi-square distribution (x ) equals zero, 

the true and estimated state vectors agree exactly. The 

2 

larger the value of x the greater the discrepancy between 
the true and estimated state vectors. 
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l xl 



What is expecced to be seen in this problem is 
a greater discrepancy in the case where the Q's of the model 
and the filter are different, due to the poor estimations of 
the state vectors. 

The above expectation seems to be the case. 

Specifically in the cases where the same Q for the filter 

and the model is used, it can be said roughly that the 
2 2 

X q - and x 375 are approximately 0.06 and 7.5 as they should 
be for two degrees of freedom (see proper tables) . 

d. Position and Velocity Error Variances (expr. b,d) 

In Cases 1 and 2 (filter has Q = 0.02) , we have 
identical position error variances and identical velocity 
error variances. After some initial variations the position 
error variance takes its steady state which is almost zero 
(expected) , while the velocity error variance takes also its 
steady state which is approximately 0.04. 

In Case 3, the increase in Q, increases P^(kjk-l)' 
and reduces P^(k|k) (see proper equations). Then the steady 
state values of P-^(k|k-l) and P^(k|k) are different and 
the result is that the position error variance oscillates 
between the values 0.0 and 0.01. For the same reason the 
velocity error oscillates between the values 0.15 and 0.24. 

The mean values of both then are greater than the corresponding 
values of Cases 1, 2. This is logical since now the Q is 
greater as it was derived! 
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Position Error Variance: 

P 11 (k|k) = P 11 (k;k-1) - G x (k) P 11 (k : k-1) (6.38) 

£ 1X { k i k— 1 ) = £]_]_ ( k— 1 | k-1) + TP 21 (k-1 | k-1) 

+ TP 12 (k-1 ) k-1 ) + T 2 P 22 (k-1 | k-1) + T 2 Q (6.39; 

Velocity Error Variance: 

P 2 2 ( k j k ) = -G 2 (k)P 2 (k|k-l) + ? 22 ( k j k- 1 ) (6.40) 

P_ 22 (k k-1) = P 22 (k-1 | k-1) + Q (6.41) 

So Q is proportional to P^ and ?22‘ 

B. AN EXAMPLE USING E.K.F. (NON-LINEAR CASE) 

1 . Problem Description 

A stationary target is located at x p = x 2 = 100- 
Bearing measurements from a moving ship are taken at 
k = 1 , . . . , N from locations as indicated. 

The following are given: 

z(k) = 0_(k) + v (k) (6.42) 

-1 X 2 

9 (k) = tan ( ' k ) k = 1 ,N (6.43) 



E [v (k) ] 



0 



(6.44) 




Movinc Shi 



Figure 6.25 Descriptive Diagram for the Non-linear Problem 



E [v (k) v ( j ) ] = o 



2 . 



'kj 



( 6 . 4 => > 



x ( 0 | 0 ) = 



10k a 2 


= 


(2 0 ) 2 


estimate 


is : 




• ^ 

x x (0 | 0) 




80 


/\ 

x 2 (0 | 0) 




120 









P (0 I 0) = 



100 



100 



N = 12 



(6.46) 



(6.47) 



(6.48) 



Based on the bearings of the moving ship, an Extended 
Kalman Filter can be used to improve the initial estimate 
of the target's position. 
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It is useful to obtain che values of x iN N ) and P(N 'A 
for observing the quality of the simulation. 

2 . Analys is 

Since the target is stationary: 



x- l (k+1) 
x 2 (k+1) 



r x (k) i 

i i 
x 2 ( k) 



(6.49, 



The model equations for non-linear problems are: 



so 



x(k+l) = a (x (k ) , u ( k ) , k ) + w(k) 



(6.50) 



z (k) = c (x (k) ) + v ( k ) 



(6.51) 



For this problem 



a (x (k) , u (k) , k) = 



(k) 
x 2 (k) 



(6.52) 



A (k) = 



3a 



dx 



(x (k | k) ) 



3x., 3x. 



I 



3x, 3x. 



3x„ 3x, 



3x, 3x, 



J 



1 0 



0 1 



(6.53) 



c(x(k)) for this problem is given by: 
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c (x (k) ) = 6 = tan 



-1 



x. 



X x - a (k) 



so 



ac 

H = d 



x ( k ) 



£0 

3x 



x(k k-1 ) 



and 



(x^ - a (k) ) 2 + x 2 
x^ - a (k) 

(x 1 - a (k) ) 2 + x 2 



(6.54) 



(6.55) 



(6.56) 



(6.57) 



To calculate the gains of the filter, the following 
equation is used: 



G (k) = P (k | k-1) H ' (k) 



I H (k) P (k | k-1 ) H ' (k) + R(k) ] 



(6.58) 



After some manipulations of the above equation in its matrix 
form, it is found that the gains are given by the following 
equations : 



P,, (k|k-l)H, (k)+P 19 (k|k-l)H~(k) 

G (k) = ±± ± if f 

H x (k) P u (k I k-1) +H 1 (k) h 2 (k) [P 21 (k I k-1) +P 12 (k I k-1) +h 2 (k) P 22 (k I k-1) +R(k) ] 

(6.59) 
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G 2 (k) — 



5 21 (k k-1 ) k, ik) ^? 22 (k , k-1; H 2 (k; 



H 1 (k) P 11 (k 1 k_1 ) (•<) h 2 (:< ) t p 2i (k I k_1) +p i2 (k k *D +K 2 (k) ? 22 (k k-l)-R(k) ] 



(6.60) 



The estimation of the states is obtained by manipu- 
lating the following equation: 

A A A 

x(k k) = x(k k-1) + G(k) [ zy k) - c (x (k k-1))] (,6.61) 



The results are: 

~ - , x (k k- 1 ) 

x-, ( k | k ) = x. (kik-1) + G, (k)z(k) - G, (k) tan * 

11 x 1 (k,k-l) - a ( k ) 

(6.62) 

A 

. x _ ( k k — 1 ) 

x 2 (k|k) = x-(k|k-l) + G~(k)z(k) - G-,(k)tan 7 

x x (k I k- 1 ) - a (k) 

(6.63) 

The error covariance matrix during the filtering phase 
of the Kalman Filter is given by the equation: 

P(kjk) = [I G (k) H (k) ]P (k | k-1) (6.64) 



By manipulating this equation as it was described previously, 
the following equations are obtained: 
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? u (kjk) = P 11 (k,k-1) [1 -G 1 (k)H 1 (k)] -P 21 (kjk-l)G 1 (k)H 2 (k) 

(6.65) 

P 12 (k k) = [1 -G 1 (k)H 1 (k)]P 12 (k k-1) -G 1 (k)H 2 (k)P 22 (k.k-l) 

( 6 . 66 ) 



P 21 (k k) = [1 -G 2 (k)H 2 (k)]P 21 (k k-1) -G 2 (k)H 1 (k)P 11 (k k-1) 

(6.67) 

P 22 (k I k) = [l -G 2 (k)H 2 (k)]P 22 (k k-1) -G 2 (k)H 1 (k)P 12 (k : k-l) 

( 6 . 68 ) 



The error covariance matrix during the prediction 
phase of the Kalman Filter is given by the equation: 

P (k | k-1 ) = A(k-l) P (k-1 | k-1) A' (k-1) +Q(k-1) (6.69) 



This yields the following results: 



P 



P 



P 



11 



12 



21 



(k | k-1) 


= P 1X (k-1 | k-1) 


(6.70) 


(k | k-1) 


= p 12 (k-1 | k-1) 


(6.71) 


(k | k-1) 


= P 21 (k-l|k-l) 


(6.72) 
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? 22 (k k-1) 



P 2 2 ( k — 1 k-1) 



( 6 . 7 3 > 



The normalized state error squared: 

[x (k) - x(k ,k) ] 'P _1 (k) (x(k) - x (k , k) ] 

is given by the same expression as in the previous problem. 

All the above results are used in a computer program 
(see Appendix A). This program basically follows the steps 
of an Extended Kalman Filter algorithm as it was described 
in Chapter III. 

From the outputs it can be seen that the position 
variance on the horizontal axis is smaller than that on the 
vertical axis. This was expected since the bearings are 
crossed with relatively small angles and in this way, there 
is a larger uncertainty on the vertical axis. 

Also from the results, it can be seen that the esti- 
mates for x^ , x 2 , are improved as time goes on and finally, 
that they are very close to the real values. 

C. MANEUVERING TARGET EXAMPLE 
1 . Problem Description 

The target is modeled by the following equation: 
x (k+1 ) = F x (k) + G w (k) (6.74) 
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This equation is discretized over time intervals of length 
T. Using Cartesian coordinates, the state is: 



and 



x = [x x y y] ' 



F = 



1 

0 

0 

0 



T 

1 

0 

0 



0 

0 

1 

0 



0 

0 

T 

1 



w = 



[w. 



w 2 ] 



G = 



T/2 

1 

0 

0 



0 

0 

T/2 

1 



(6.75) 



(6.76) 



(6.77) 



(6.78) 



Eiw(k) ] 



0; E[w(k) w * ( j ) J 



Q6 R . (6.79) 



The initial estimate is x(0[0) with covariance P(0|0). 
It is assumed that only position measurements are 
available : 



z(k) = Cx(k) +v(k) 



(6.80) 



where 
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c = 



10 0 0 
0 0-10 



(6.81) 



E [ v (k) ] = 0; E [v (k) v' ( j 



= Ro, . 

- k] 



6.32; 



The following are given: 



T = 10s, Q = 0, R n - R 22 = 10 4 m 2 , and 



R-j^ 2 = 500m 



The initial conditions of the target are: 



x ( 0 ) = 200m, x ( 0) = 0, y(0) = 10,000m, y(0) = -15m/s 



The target is on a constant course and speed until 

time t = 400s, when it maneuvers a slow 90° turn in the x 

xv 2 

direction with acceleration inputs u = u J = 0.075m/s . It 

completes the turn at t = 600s (after 20 sampling times) and 

from then on the accelerations are zero. The second turn, 

also of 90°, is fast. It starts at t = 610s with acceleration 
2 

of 0.3m/s and is completed at t = 660s, i.e., after 5 
sampling times. 

A simulation of the maneuvering target can be done 
in x coordinate only (the results for the y coordinate are 
similar) , using the method of "White noise model with adjusta- 
ble level . " 
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The quality of the performance of the target's simu- 
lation can be observed from the following plots: 

a. True and estimated trajectoris of the target. 

b. Monte Carlo runs of position root mean square 
error (r . m. s . e. ) : 




50 

L 

i=l 




(k | k) ] 2 



1/2 



c. The same as b. for the velocity r.m.s. error. 

2 . Analysis 

The target's movement is presented in Figure 6.26. 




Figure 6.26 Movement of the Manuevering Target 
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For treating this problem, a linear Kalman Filter 
can be used. The only difference is that the actual trajec- 



of one. The different phases of movement, of the target A,B, 
C,D 7 E, can be formed in proper mathematical formulas, which 
can feed the position measurements one after the other, 
depending on the sampling time. 

The initial estimates of the states in the x direction 
were obtained as : 



tory of the target is now given by several formulas instead 



x (0 . 0) 



x L (0 .0) 




x ( 0 | 0 ) = x 2 ( 0 , 0 ) = [z 1 (0) - z ] _ (-1) ]/T 



where 




x ( 0 ) - x ( 0 ) T + v 1 (-l) 



(6.35) 




x(0) + v^ ( 0 ) 



( 6 . 86 ) 



v ~ N(0,R u ) 



(6.87) 



The initial covariance matrix is then: 



R 



11 




P(0 | 0) 



( 6 . 88 ) 



R ll /T 2R 11 /T 



2 
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The rest of the calculations needed for developing 
the Kalman Filter algorithm are very similar to those 
developed in the first example. 

If the maneuvers are not taken into account and the 
filter runs as usual, the tracking results are very poor 
during the maneuvers (see pp . 114-115). 

Using the "White noise model with adjustable level" 
method, a self-ad j ustment of the filter can be accomplished 
as follows: 

After the square norm of innovations: 

£ v (k) = V ' (k) S _1 v (k) 

exceeds a chosen value, the target is considered as maneuver- 
ing and a change (increase) in Q, is made. The increase 
in Q increases the gain G, so the state estimation: 

x (k | k) = x (k | k-1) + G(k) [z (k) - c(k)x(k|k-l) ] (6.39) 

is more dependent on the second term since the state estimation 
x(k|k-l) is less reliable during the maneuver. The filter 
continues its run with the new value of Q' until has again 
a lower value than that of the chosen threshold. Then the 
filter runs again with the given initial value of Q (Q = 0.0) . 

After the addition of the above provision in the 
tracking computer program, the results were much better than 
those of the non-ad j ustable filter. 
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Three categories of results were obtained. 



In each 



one, different values of Q's and threshold levels were used. 
The three categories are: 

a. Observed and estimated trajectories versus sampling 
time (50 Monte Carlo runs). 

b. Observed and estimated trajectories versus sampling 
time (single runs) . 

c. Mean square position error (50 Monte Carlo runs . 

The above results (plots) can be seen on pp . 116-149 
A tabulation and qualification of these results has been 
developed as follows: 



Symbols 




S : 


Constant course and speed 


1M: 


First Maneuver 


2M: 


Second maneuver 


Very Good: 


The estimated trajectory follows the time trajector 
with great reliability (coincidence of the paths). 


Good : 


The estimated trajectory is still very reliable, 
i.e., follows the true path closely but they 
don't coincide exactly. 


F. Noise: 


The estimated trajectory follows the noise. 


F . Noise E . : 


The estimated trajectory follows the noise 
exactly . 


Deviation : 


There is deviation between the true and estimated 
tra j ectory . 


S . Deviation : 


There is strong deviation between the true and 
estimated trajectory. 
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TABLE II 



Observed and Estimated Trajectories 
(50 Monte Carlo Runs) 



THR. 

Q 

S 

0.1 1M 

2M 

S 

3 1M 

2M 

S 

10 1M 

2M 

S 

100 1M 

2M 

S 

1000 1M 

2M 



0.1 

Very Good 
Deviation 
S. Deviation 

Very Good 
Good 

Deviation 

Very Good 
Very Good 
Good 

Very Good 
Very Good 
Good 

Very Good 
Very Good 
Very Good 



3 

Very Good 
S. Deviation 
S. Deviation 

Very Good 
Deviation 
Deviation 

Very Good 
Good 

Deviation 

Very Good 

Good 

Good 

Very Good 
Very Good 
Good 



20 

Very Good 
S. Deviation 
S. Deviation 

Good 

Deviation 

Deviation 

Good 

Deviation 

Deviation 

Good 

Deviation 

Deviation 

Deviation 
Deviation--F . 
Deviation — F . 



Noise 

Noise 
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TABLE III 



THR. 

Q 



0 . 1 



3 



10 



100 



1000 



Observed and Estimated Trajectoris 
(Single Runs) 





i — 1 

o 


3 


20 




s 


Very Good 


Very Good 


Very Good 




1M 


Deviation 


S. Deviation 


S. Deviation 




2M 


S. Deviation 


S. Deviation 


S. Deviation 




S 


F. Noise 


Very Good 


Very Good 




1M 


Good 


Good 


Deviation 




2M 


Good 


Deviation 


Deviation 




S 


F. Noise E. 


F. Noise 


Good 




1M 


Good 


Good 


Deviation 




2M 


Good 


Good 


Deviation 




S 


F. Noise E. 


F. Noise 


Good 




1M 


F. Noise E_. 


Good 


Deviation- -F . 


No 


2M 


F. Noise E. 


F. Noise 


De viation--F . 


No 


S 


F. Noise E. 


F. Noise 


F. Noise 




1M 


F. Noise E. 


Good 


Deviation--F . 


No 


2M 


F. Noise E. 


F. Noise 


Devia t ion--F . 


No 
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Symbols 



1MAX : 


Maximum value of m.s.e. during the first 
maneuver 


2 MAX : 


Maximum value of m.s.e. during the second 
maneuver 


3MIN: 


Minimi um value of m.s.e. during movement 
under constant course and speed. 


MAN: 


Maneuvers 
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Root Mean Square Position nrror (X Direction-) 



(50 Monte Carlo Runs) 



THR. 

Q 




i — 1 

o 


3 


20 




1MAX 


227 


370 


730 


0 . 1 


2 MAX 


257 


300 


370 




3MIN 


32 


30 


30 




1MAX 


82 


147 


4 5 5 


4 


2 MAX 


108 


190 


365 




3MIN 


48 


37 


30 




1MAX 


84 


130 


425 


10 


2 MAX 


92 


168 


340 




3MIN 


52 


37 


40 




1MAX 


102 


107 


375 


100 


2 MAX 


98 


135 


305 




3 MIX 


6 4 


46 


30 




lMAX 


114 


102 


345 


1000 


2 MAX 


110 


122 


295 




3MIN 


72 


51 


30 




1MAX 




102 




3000 


2 MAX 
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3MIN 




52 





TABLE V 



Conclusions From Tables II, III, IV 



TABLE II 



TABLE III 



High Q 



S Good 



MAN Very Good 



F. Noise 
F. Noise E. 
F. Noise 
F. Noise E. 



Low Q 



S Very Good 



MAN S. Deviation 



Very Good 

Deviation 
S. Deviation 



S Deviation 

High 

Threshold 

MAN Deviation 



Very Good 
Good 

Deviation 
S. Deviation 



S Very Good 

Low 

Threshold 

MAN Very Good 



F. Noise 
F. Noise E. 
Good 

F. Noise 



TABLE IV 

r .m. s . e . 
increases 
r . m. s . e . 
decreases 

r .m. s . e . 
decreases 
r . m . s . e . 
increases 

r . m. s . e . 
decreases 
r . m . s . e . 
increases 

r . m. s . e . 

increases 
r . m. s . e . 
decreases 
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50 Monte Carlo Huns 



08S. RND ESTIR.TRRj. NON -RD JUST. FILTER 

SINGLE RUN OF THE PROGRRtl 
CONT1N. LINE: MEASUREMENTS TRRjECT. 

DRSHCD LINE: F ILTCR S ESTIM. TRfUECT . 
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Figure 6.28 Obs. and Fstim. Traj. Non 
' Single Run of the Program 
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50 Monte Carlo Runs 
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Fiijure 6.30 Observed and EsLim. Trajecl. O' 100 TK 
50 Monte Carlo Runs 
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ligure 6.31 Observed and bsfcim. Traject. O' 
50 Monte Carlo Runs 



OBSERVED PND EST1M. 7RPUECT. C-.l ' 

SINGLE RUN OT THE PROGRAM 
CONTiN. LINE: MEASUREMENTS TRAjECT. 
DASHED LINE: E 1LTER S ESTiH. TRBjECT. 
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Figure 6.32 Observed and Estim. Tra jee 
50 Monte Carlo Runs 



OBSERVED RND ESTih. TR8JLCT. 0-3 TR- 

SINGLC RUN OT THC PROGRPfl 
CONT'N. LINE: NCPSURCMCNT5 TRfUCCT. 
DRSMf.0 LINC: riLTER S lSTIM. TRfljCCT . 







1.5'jL.J.r “C_ n.I'J '_'Nd jDAireSuC 
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figure 6.33 Observed and list ini. Ttaje 
Single Run of the Progiam 



OBSERVED PND EST1M. 7RPJECT. 0-10 TR- 
S INGLE RON or THE PROGRPH 
CONTiN. l!NE: MEASUREMENTS TRRjECT . 
DRShED LINE: DETER 5 ESTIM. TRRjECT. 
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F iyure 6.34 



QUSCRVCD AND CST1M. TRPJECT . 0-100 TR- 

S INGLE RUN Or TML PROGRAM 
00NT1N. LINE: NEPSUREMEWS 7RPJEC7. 
DHSHEC LINE: flLTER S C57IM. 7RRJEC7 . 
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Observed and h'st. ini. 



OBSERVED PND ESTIM. TRPJECT. C-IOOO TR 

SINGLE RUN Of THE PROGRAM 
CONTIN. LINE: MEASUREMENTS TRRjECT . 
ORSHEO LINE: FILTER 5 E5TIN. TRAjECT . 
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Figure 6.36 Observed and Est.Lm. Trajeet. O' 1000 TR 
Single Run of the Program 



OBSERVED FIND CSTIh. TRflJEXT. O-.l TR-3 

SINGLE RUN or THE PROGRAM 
CONTIN. LINE: MEASUREMENTS TRAjECT. 

OASHEO LINE: MLTER S EST1M. TRAJCCT . 
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rigure 6.37 Observed and Rstim. Tin jo 
Single Run of the Program 



OBSERVED AMD ESTlh. TRRJECT . 0-3 TR-3 

SINGLE RUN Of THC PROGRAM 
CONTIN. LINE: MEASUREMENTS TRflJECT. 
DASHED LINE: TIL TER S ESTIM. TRflJECT . 
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Figure 6.38 Observed and CstLin. Traject. 

Single Run of the Program 



OBSERVED AND ESI I A. 7RAJECT. O-IO TR-3 

SINGLE RUN Gf THE PROGRAM 
CONTIN. LINE: MERSUREPENTS TRflaECT . 

ORSHEO LINE: TILTER S ESTIP. TRfiJLCT. 
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I icjure 6.39 Observed and tistim. 



OBSERVED PND ESTIM. TRRJECT. O-IGO TR-3 

5 INGLE RUN OT THE PROGRAM 
CONTIN. LiNE: MEASUREMENTS TRAJECT. 

OASHEO LINE: KILTER S ESTIM. TRfUECT. 
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Figure 6.40 Observed and Estim. Traject. 

Single Run of the Urogram 



OBSERVED HMD EST1K. TRBJECT . 0-1 GOO TR 

SINGLE RUN or THE PROGRAM 
CONTiN. LINE: MlOSUREMENTS TRfiJECT . 
DfiSHED LINE: TILTER S EST1M. TRRjECT. 
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Figure 6.41 Observed and Estim. Trajecl . O' l 00() Tl< 
Single Run of the Program 



OBSERVED RN9 EST1M. TRflJECT. 0- . I ' 

SINGLE RUN Of THE PROGRAM 
CONT ! N . LINE: MEASUREMENTS TRAJECT. 
DASHED LINE i FILTER 5 ESTIM. TRAjECT. 
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Figure 6.4 2 Observed and Estini. 



OBSERVED AND ESTlh. TRflJECT. C-3 TR-20 
SlNGt.C RUN or THt PROGRAM 
CONTiN. L1NC: NCOSuRl MCNTS TRRjfXT. 

ORSlCO UNt': flLTtR S PSTIM. TRfl.JCCT . 




130 



icjure 6.4 3 Observed and Kst ini. 



OBSERVED RND EST1M. IRflJECT . 0- 1C TR-20 

SINGLE RUN OT THE PROGRAM 
CONTIN. LINE: MEASUREMENTS TRAjECT. 

DASHED LINE: TILTCR S CST1M. TRAjECT. 
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Figure 6.44 Observed and Estim. Trajocl . Q’ 10 TR 
Single Run of the Program 



OBSERVED PND ESTIM. TRPJECT . O-lOO TR-20 
SINGLC RUN or THE PROGRRN 
CONT!N. LINE: MEASUREMENTS TRPJECT. 

OPSHED LINE: flLTER S EST !M. TRPJECT. 
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ure 6.45 Observed and Estim. Traje 
Single Run of the Program 



OBSERVED AND EST1R. TRRJECT. O-IOCO TR-20 

SINGLC RUN Of THE PROGRAM 
CONTIN. LINE: MEASUREMENTS TRAJECT . 

OASHEO LINE: FILTER S ESTiM. TRAJECT. 
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Figure 6.46 Observed and Estini. Traject 
Single Run of the Program 
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Figure 6.47 Mean Square Frroi of Pos. 



MEAN SQUARE ERROR OF POS. 
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Figure 6.4 8 Mean Square Error ol I’os. 



ICfiN SOURRE ERROR OR POS. 
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Figure 6.4 9 Mean Square Error of I’os . O' 10 TK 



MEAN SQUARE ERROR OF POS. 0= 1 00 TR- 
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Figure 6.50 Mean Square Error of kos . O' 100 Tk 
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Figure 6.51 Mean Square Error of Pos . Q* - 1000 Tk 



MEAN SQUARE ERROR OF POSIT. 
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Figure 6.52 Mean Square Error of Pos. 



MEAN SQUARE ERROR OF POSIT. 0-4 TR=3 
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Kicjure 6.53 Mean Squar e ttrror of l‘os . 



MERN SOURRE ERROR OF POSIT. 0=10 TR= 
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Figure 6.54 Mean Square Error of Pos. 
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SAMPLING TIME 



MEAN SOURRE ERROR OE POSIT. 0=1000 TR=3 
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Figure 6.56 Mean Square Error of Fo s. Q ' 1000 Tl< 



MEAN SQUARE ERROR OF POSIT. 0=300 
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Figure 6.57 Mean Square Error of Fos, Q' = iOOO TU 



MERN SOURRE ERROR OF P05. 0- . i TR = 20 
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Ficjure 6.58 Mean Square Error of Fos. 
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•'icjure 6.59 Mean Square 
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figure 6.60 Mean Square Error of Eos. 



MEAN SQUARE ERROR OF P05. 0=100 TR-20 
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Figure 6.61 Mean Sijuare Frroi of Fos. 
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Figure 6.62 Mean Square Error of Pos. Q 1 1000 TU 



Since the observed and estimated paths in the Monte 
Carlo case are a-veraged, the tracking results are better 
(smoother) and it is not easy to observe the various affec- 
tions of extreme values of Q and the threshold level. For 
this reason, only three characteris tic plots of this case 
are presented and the following conclusions are based mainly 
on the other two cases: 

a. When the target is moving with constant course and 
speed, it is desirable to have low values of Q and 
high threshold level. 

b. 'When the target is maneuvering, it is desirable to 
have medium or high values of Q and low threshold 
level . 

Based on the above results, an attempt was made to 
run the filter using three levels of process noise Q (two 
threshold levels) . In this way the two types of maneuvers 
(slow-fast) were treated separately. The threshold levels 
chosen were 1 and 5. The mean square error of position and 
velocity in the x direction over 50 Monte Carlo runs was 
plotted. The results (pp. 152-171) were very good especially 
for the position m.s.e. A tabulation of these results follows 
as Table VI. 
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TABLE VI 



Multiple 



THR. 

Q / 



1MAX 

0,0.1,100 2 MAX 

3MIN 

1MAX 

0,3,100 2 MAX 

3MIN 

1MAX 

0,10,100 2 MAX 

3 MIN 

lMAX 

0,0.1,10 2 MAX 

3MIN 

lMAX 

0,1,10 2 MAX 

3 MIN 

lMAX 

0,3,10 2 MAX 

3MIN 

lMAX 

0,0.1,1000 2 MAX 

3MIN 

lMAX 

0,1,1000 2 MAX 

3MIN 

lMAX 

0,3,1000 2 MAX 

3MIN 

lMAX 

0,10,1000 2 MAX 

3MIN 



Process Noise Levels 



R.M.S. Position Error 
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R.M.S. Velocity Error 
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sampling time 



MEAN SOUARE POSITION ERROR 
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Figure 6.64 Mean Square PosiLion Error 
Q ' Levels : 0,3,100 
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Mean Square Posit ion 



MEAN SOURRE POSITION ERROR 
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SAMPLING TIME 
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Figure 6.61 Mean Square Position Frror 
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Figure 6.68 Mean Square PosiLion 
Q ' Levels : 0 , 3 , i 0 
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SAMPLING TIME 



MERN SOUflRE POSITION ERROR 
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SAMPLING TIME 
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SAMPLING TIME 
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SAMPLING TIME 

Figure 6.72 Mean Square Position Frr 
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SAMPLING TIME 



MEAN S0UARE VELOCITY ERROR 
MANY 0 LEVELS 

0 LEVELS: O'3-IOO 
TORES. LEVELS: l~S 




aoyyo Q33dS swy 



163 



Mean Square Velocity 
Q* Levels: 0,3,100 
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SOMPLING TIME 
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figure 6.76 Mean Square Velocity fr 
O' Levels: 0 , . 1 , 1 0 
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'lyure 6.77 Mean Square Velocity 
O' Levels : 0,1,10 
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HERN SOURRE VELOCITY ERROR 
MANY 0 LEVELS 
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Figure 6.79 Mean Square Velocity Frroi ; Many (.) ‘ 
Q' Levels : 0,1,1000 
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SAMPLING TIME 
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SRRPLING TIRE 



MEfiN SOUflRE VELOCITY ERROR 
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Figure 6,82 Mean Square Velocity 



3 . Comparison With Other Approaches 



Bar-Shalom [Ref. 9], has developed a different approac 
for the same problem. In this case the tracking filter 
operates in its normal mode in the absence of any maneuver. 
Once a maneuver is detected, a different stare model is used 
by the filter. The acceleration is added as a new state 
component. The extent of the maneuver as detected, is then 
used to yield an estimate for the extra state component, 
and corrections are made on the other state components. The 
tracking is then done with the augmented state model until 
it will be converted to the normal model by another decision. 

The results of this approach can be seen in Figure 
6.83. The curve having the indication "VD" belongs to the 
above algorithm while the other curve (indication "IE") 
belongs to a different approach developed by Chan [Ref. 10] . 
The outline of the latter algorithm is the following: 

When a maneuver is detected, the magnitude of the 
acceleration is identified in a least squares format. The 
result is used in conjunction with a standard Kalman Filter 
to estimate the state of the vehicle. The aim of the accelera 
tion input estimation is to remove the filter bias caused 
by the target deviating from the assumed constant velocity, 
straight line motion. 

It can be easily see (Figure 6.83) that the algorithm 
of Reference 9 is superior to that of Reference 10. 

The "white noise model with adjustable level" method, 
modified with multiple levels of process noise Q (one for 
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Figure 6 . B 3 References 



each type of maneuver) , can be compared with the algorithms 
of References 9 and 10. Specifically this paper's algorithm 
is better than both the algorithms in the r.m.s. position 
error results and worse in the r.m.s. velocity error results. 
In the r.m.s. velocity error case, the approximate maximum 
value of the two other algorithms is 11, while this paper's 
algorithm maximum value is 13.3 (not too big a difference' . 

In the r.m.s. position error case, the maximum values for 
References 9 and 10 algorithms are 200 and 125, respectively, 
while this paper's algorithm has achieved an approximate 
value of 95. 
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VII. CONCLUSIONS 



In the examples presented, it was seen that the Kalman 
Filter has the ability and flexibility to treat various 
cases of the target tracking problem with very sacisf accory 
results . 

The Kalman Filter is always initialized by the user, 
providing the initial estimate x(k) and its corresponding 
estimate error covariance matrix P(k). The initialization 
is of great importance for the filter's performance. A poor 
initialization needs more observations and time for the 
algorithm estimate to converge toward the value of the state 
vector . 

In all cases the importance of the gain matrix G(k) of 
the Kalman Filter was significant, i.e., since it is inversel 
proportional to time, it weights the correction term 
[z(k) -x(kjk-l)] less heavily as time progresses, and so the 

A 

state estimation x(k|k) depends more on the state estimation 
of past time x(k|k-l). 

It was also noted that the Q matrix accounts for any 
model inaccuracies. For a filter in a steady state condition 
the Q also serves to prevent the gain matrix G(k) from 
approaching zero by always insuring uncertainty in the pre- 
dicted error covariance matrix. Q is also the key variable 
for treating the maneuvering target problem with the "White 
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noise model with adjustable level" method. As it was ex- 
plained in the corresponding chapter, variations in Q cause 
proportional variations to the gain G(k) , which weights the 
correction factor in the equation: 

x ( k : k ) = x(k k-l) M G(k)[z(k) - c ( k ) x ( k k- 1 ) ] 

Larger values of Q tend to cause the state estimate to favor 
the observation. 

The error covariance matrix and specifically some of 
its elements (position and velocity error variance) , were 
very useful as indicators of the filter's performance. 

In the case of the maneuvering target, the "White noise 
model with adjustable level" method was very reliable and 
comparable to the algorithms of References 9 and 10. Addi- 
tionally, it is much simpler than the other two methods. The 
only disadvantage was the greater value of the velocity r.m.s. 
error. In the example presented, the difference was not 
significant but it is estimated that it would increase in 
cases where the target maneuvers more drastically (greater 
accelerations) . 

The Monte Carlo simulation should be used to provide 
statistical results for any stocahstic process which is 
represented by pseudorandom numbers. This simulation is a 
rigorous statistical procedure, to compare two or more 
algorithms . 
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More specific conclusions concerning the examples presented 
have already been provided in the corresponding analyses of 
them. 

It is suggested that a continuation of this research 
should include the following: 

a. Problem of maneuvering target which uses greater 
acceleration driving the maneuvers. The "White noise 
model with adjustable level" method should be applied and 
the results in position and velocity r.m.s. errors should be 
obta ined . 

b. Problem of multiple targets under various situations 
(maneuvering, in clutter, etc.). 

c. Problem of single and multiple maneuvering targets 
in clutter. It is an interesting and important topic not 
included in this paper due to lack of time. 

d. Study the implementation of a Kalman Filter algorithm 
using special purpose hardware. 
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